/////////////////////////////////////////////////////////////////////
// \author Jeongnim Kim
//
// An example to use TinyVector and Vector classes
// To use SPRNG library:   "make sprng.x"
// To use C rand function: "make random.x"
// Testing generation of TinyVector<T,3> and vector<TinyVector<T,3> >
//////////////////////////////////////////////////////////////////////
#include <iostream>
#include <vector>
#include <iomanip>
using namespace std;

#define DEBUGMSG(msg)

#include "Utilities/RandomGenerator.h"
#include "OhmmsPETE/TinyVector.h"
#include "ParticleBase/ParticleAttrib.h"
#include "ParticleBase/ParticleUtility.h"
#include "Utilities/Timer.h"

int main(int argc, char **argv) {

  typedef TinyVector<double,OHMMS_DIM> Pos_t;
  typedef ParticleAttrib<Pos_t> PosArray_t;

  bool makeseed = true;
  Random.init(makeseed); 
  RandomVector<Pos_t, RandomGenerator_t> rvgen(Random);

  int n=10;
  int i=0;
  while(i<argc) {
    if(!strcmp(argv[i],"-n")) {
      n = atoi(argv[++i]);
      break;
    }
    i++;
  }

  cout << endl << "Starting vector operations.\n";
  Timer myclock;
  myclock.start();  
  PosArray_t R(n), P(n), Q(n);
  myclock.stop();
  cout << "Creating vectors = " << myclock.cpu_time() << endl;

  myclock.clear();
  myclock.start();
  RandomSequence<PosArray_t, RandomGenerator_t>::apply(R, rvgen);
  RandomSequence<PosArray_t, RandomGenerator_t>::apply(P, rvgen);
  myclock.stop();
  cout << "Random vector assignment = " << myclock.cpu_time() << endl;

  myclock.clear();
  myclock.start();
  Q = 1.5*P + 2.1*R + 1.8*R;
  myclock.stop();
  cout << "Operation by expression template= " << myclock.cpu_time() << endl;
  myclock.clear();
  myclock.start();
  for(int i=0; i<n; i++)
    Q[i] = 1.0*P[i] + 2*R[i];
  myclock.stop();
  cout << "Operation by a for-loop = " << myclock.cpu_time() << endl;
 
  if(n<=100) {
    for(i=0; i<n; i++)   {
      cout << R[i] << " " << P[i] << " " << Q[i] <<  endl;
    }
  }

  typedef ParticleAttrib<double> ValArray_t;
  ValArray_t a(n), b(n), c(n);
  for(int i=0; i<n; i++) b[i] = R[i][0];
  for(int i=0; i<n; i++) c[i] = P[i][0];

  cout << endl << "Starting scalar operations.\n";
  myclock.clear();
  myclock.start();
  a = 1.5*b + 2*c;
  myclock.stop();
  cout << "Operation by expression template = " << myclock.cpu_time() << endl;

  myclock.clear();
  myclock.start();
  for(int i=0; i<n; i++)
    a[i] = 1.0*b[i] + 2*c[i];
  myclock.stop();
  cout << "Operation by a for-loop = " << myclock.cpu_time() << endl;
  
}

/***************************************************************************
 * $RCSfile: vector.v0.cpp,v $   $Author: jnkim $
 * $Revision: 1.1.1.1 $   $Date: 2004/07/23 21:04:02 $
 * $Id: vector.v0.cpp,v 1.1.1.1 2004/07/23 21:04:02 jnkim Exp $
 ************************************************************************/
